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ABSTRACT 

Using a new parallel algorithm implemented within the Visit framework, we analysed 
large cosmological grid simulations to study the properties of baryons in filaments. The pro¬ 
cedure allows us to build large catalogues with up to ~ 3 • 10 4 filaments per simulated volume 
and to investigate the properties of cosmic filaments for very large volumes at high resolution 
(up to 300 3 Mpc 3 simulated with 2048 3 cells). We determined scaling relations for the mass, 
volume, length and temperature of filaments and compared them to those of galaxy clusters. 
The longest filaments have a total length of about 200 Mpc with a mass of several 10 15 Mq. 
We also investigated the effects of different gas physics. Radiative cooling significantly mod¬ 
ifies the thermal properties of the warm-hot-intergalactic medium of filaments, mainly by 
lowering their mean temperature via line cooling. On the other hand, powerful feedback from 
active galactic nuclei in surrounding halos can heat up the gas in filaments. The impact of 
shock-accelerated cosmic rays from diffusive shock acceleration on filaments is small and the 
ratio of between cosmic ray and gas pressure within filaments is of the order of ~ 10 — 20 
percent. 

Key words: galaxy: clusters, general - methods: numerical - intergalactic medium - large- 
scale structure of Universe 


1 INTRODUCTION 


Simulations of structure formation predict that about 
50% of the baryons reside in a filamentary web of ten¬ 
uous matter at temperatures f 1 0 5 — 10 7 K connecting 


already virialized structures (e.g. Cen & Ostriker 1999 
Dave et al. 20011. This plasma is referred to as ‘warm 


Hot Intergalactic Medium” (WHIM). 

The total matter (mostly Dark Matter, DM) of the 
cosmic web is traced by galaxies that populate filamen- 


sizes in excess of 100 Mpc/h (e.g. Zeldovich et al. 

1982 

Einasto et al. 1984; Geller & Huchrall989; Gott 

et al. 

2005 

). A first glimpse of this cosmic web in the 


slice (e.g. de Lapparent et al. 19861). In recent years, 


this view has been expanded by the 2dFGRS (e.g.|Col 


less et al.|2003[), SPSS (e.g . Tegma rk~et al.||2004| i and 


IRAS (e.g. |Courtois et al.| 2013) galaxy redshitt sur¬ 

veys. 


A direct detection of the gas in the cosmic web is 
more challenging because the low densities and tem¬ 
peratures are unfavourable, and so far only very little 


evidence has been produced. Baryons in filaments may 
be revealed in the soft X-ray band but the few reported 
detections are still controversial (e.g.|Finoguenov et al. 


2003 

Werner et al. 2008 

Nicastro et al. 2010,12013). 

In the 
have b 

radio band, only very few possible detections 

een reported (Bagchi et al. 2002; Kronberg et al. 

2007 

Giovannini et al. 2010, 

Farnsworth et al.,2013). 


These sources are more likely to be related to merger 
shocks than to the accretion shocks described by nu¬ 
merical simulation. 


Recently, the Sunyaev-Zeldovich (SZ) effect has 
been used to probe filamentary gas connections be¬ 
tween galaxy clusters, and a first indication for the de¬ 
tection of a filament between the cluster pair A399- 


A401 has been reported by Planck Collaboration et al. 
([2013J). They observed a significant thermal SZ signal 
in the regions beyond the virial radii. A joint X-ray SZ 
analysis constrained the temperature of the filamentary 
region to kT = 7.1 ± 0.9 keV and the baryon density to 
3.7 ± 0.2 • 10 -4 cm -3 . However, this should represent 
only the dense and small-size version of much larger 
(~ 10 — 100 Mpc long) objects that the cosmic volume 
can contain. 
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The evolution of the cosmic web can be explained 
through the interaction of the initial pattern of density 
waves at different scales, with random and uncorre¬ 
lated spatial phases. The main skeleton of the cosmic 



information embedded in the pattern of initial cosmo¬ 
logical fluctuations. During their evolution initial fluc¬ 
tuations interact in a non-linear way, with the genera¬ 
tion of non-random and correlated phases, which lead 


to th e sp atial pattern of the present cosmic web (Bond 
et al .|1996 1. The non-linear evolution can be described 
by N-body simulations, that provide an accurate pic¬ 
ture of the evolution of the gravitational potential and 
of the Dark Matter halos. 

With the growing capabilities of N-body cosmo¬ 
logical simulations, increasing effort has been devoted 
to the implementation of reliable methodologies for 
the identification of complex structures in the matter 
distribution, aiming at the accurate segmentation of the 
cosm ic web into cluster s, filaments, walls and voids 
(e.g. |Cautun et al.||2014| and references therein for a 
recent review). The list of attempted methods is too 
long to enable a complete summary here. The algo¬ 
rithms can be broadly grouped into: a) geometrical and 
tessellation methods, based on the topological analysis 


ematical approaches (e.g. Stoica et al. 2005; Sousbie 

et al. 2008; 

Aragon-Calvo et al. 

2010; Gonzalez & 

Padil a 2010 

Bleuler et al. 2014 

Chen et al. 2015), 


b) morphological methods, that classify the 3D distri- 
bution of matter based on the density Hessian and/or 


the tidal or velocity shear fields (e.g. Hahn et ah 
2007} Aragon-Calvo et aT7||2QQ7t |Cautun et al.||2U 1 4^7 


rag< 

Recently, |Cautun et al.|(|2014|) compared the results of 
different filament detection methods with the NEXUS 
algorithm, which is an advanced multiscale analysis 
tool to identify the different morphologies of cosmic 
structures, based on the analysis of density, tidal field, 
velocity divergence and velocity shear as tracers of 
the Cosmic Web. They concluded that most methods 
agree on the largest filaments, but parameters such as 
the diameter of filaments depend on the resolution and 
the method. Most of the differences are found in the 
most rarefied environments, where density contrasts 
are very small and methods working on shears or local 
differences are less sensitive. As a result, different 
algorithms agree well on the total mass in filaments 
(~ 90 percent of the total mass is captured in all 
methods) and worse on the total volume (only a ~ 60 
percent of the volume is equally identified). 


So far, the study of the properties of the baryons 
in the low-density components of the cosmic web 
(such as in filaments) have been investigated less 
systematically. Low-density environments outside of 


galaxy clusters have been studied with cosmologica 

hydrodynamic al simulations (e.g. 

Dolag et al. 

2008 

Dave et al. 

2001, 

Viel et al. 

2005 Dolag et al. 

2006 


and references therein). The effect of energy feedback 


from galactic activity on the obs ervable propertie s of 
the WHIM has been an alysed by Cen & Ostriker (e.g 


2006) ; IRoncarelli et al.| ( e.g |2 006); Kang et aid (e.g 

2007) ; Hallman et al. (e.g 2007); Roncarelli et al. (e.g 


2012); Bolton et ai.|(e.g|2014|);|Vogelsberger et ai.|(e.g 
2014). Hydrodynamic al simulations suggest that re¬ 


gions of moderate overdensity such as filaments, host a 


significa nt fraction of th e WHIM (e.g. Cen & Ostriker 


1999 Dave et ah] 2001 1 , and they are surrounded by 


strong stationary accretion shocks, where the cosmic 


Miniati et al. 

2000; Ryu et al. 

2003; 

Pfrommer et al. 


(2012) simulated the impact of different physical 
processes as well as of the scale dependencies in 
the formation of ~ 5 Mpc filaments starting from 
a single-scale perturbation. The simulated filaments 
exhibit an isothermal core, whose temperature is 
balanced by radiative cooling and heating due to the 
UV background, yielding a multiphase medium. The 
WHIM is als o expected to host supersonic turb ulent 
motions (Iapichino et al.|2011| Vazza et al.|2014a|), and 
the decay of such motions can cause the amplification 
of we ak primordial magnetic fields up to the ~ nG 
level ( Ryu et al.|2008a ; Vazza et al.|2014a I. 


In this work, we combined state-of-the art cosmo¬ 
logical hydrodynamic simulations performed using the 
ENZO code (Section® with a novel methodology for 
the identification of the cosmic web, based on the gas 
matter distribution. The gas distribution can be accu¬ 
rately described by ENZO ’s Eulerian hydrodynamic 
solver irrespectively of the mass density. Using the gas 
mass density, we have developed an Isovolume based 
technique to identify the cosmic web structures. Iso¬ 
volumes represent a class of algorithms that tackle the 
problem of separating regions of space with distinct 
physical properties, in our case with density above and 
below a given threshold. We tuned the selection in gas 
overdensity specifically to extract the mass distribution 
associated with the filamentary structure of the cosmic 
web. However, the same approach can also be used to 
identify other structures, such as voids, sheets and ha¬ 
los. 

The resulting methodology is described in detail in 
Section [3] The corresponding performance and scaling 
properties on high-performance computing systems are 
presented in Appendix A, together with a discussion on 
the influence on the results of the parameters character¬ 
izing our approach. 

The methodology is simpler than the aforemen¬ 
tioned methods while avoiding the drawbacks of al¬ 
gorithms that rely on DM particles. In regions where 
the matter density is comparable to the mean cosmic 
density, as in the WHIM, the numerical noise arising 
from the graininess of the DM particles distribution 
can affect the estimates of the properties of the cos¬ 
mic web. Further advantages of using the gas compo¬ 
nent are that it is immediately related to observations, 
and that the effect of additional physical processes (e.g. 
cooling and feedback from galactic activity) influenc¬ 
ing the gas properties, can be studied in detail. Using 
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the gas component as a tracer of the gravitational po¬ 
tential, the method is insensitive to cosmic shear fields 
and to the local dynamics of DM. 

We have analysed the statistical and thermody¬ 
namic properties of filaments identified through our 
procedure, and investigated their dependence on the 
spatial resolution and the assumed physics, in partic¬ 
ular the cooling and hea ting processes. The results 
are presented in Sections |3.3| and [4] and discussed in 
Section [5] Our finally summary is given in Section [6] 
where we discuss possible future numerical develop¬ 
ments that will allow a more quantitative prediction of 
the observable properties of the WHIM through the in¬ 
clusion of mechanisms such as chemical enrichment, 
metal-dependent cooling and magnetic fields. 


2 SIMULATIONS 

Our simulations are run with a customised version of 


in 


th e grid code ENZO (Bryan et al. 2014), presented 


Vazza et al. (201 4b|) . ENZO is a parallel code for 


cosmological (magneto-)hydro-dynamics, which uses 
a particle-mesh N-body method (PM) to follow the 
dynamics of the DM and (in this case)_ the Piecewise 
Parabolic Method (PPM, |Colella & Woodward[[1984] > 
to evolve the gas component. On the basis of the pub¬ 
lic version of ENZO , we have implemented specific 
methods to describe the evolution and feedback of 
cosmic-ray (CR) particles (Vazza et aLl|2012|, as well 
as our implementa tion of energy release from AGN 
(Vazza et ah|2013j) and supernovae. The suite of sim¬ 
ulations analysed in this work has been designed to 
study the properties of CRs, the acceleration of rela¬ 
tivistic hadrons at shocks and their energy feedback on 
the baryon gas (Vazza et al.|2014bj). 

The simulations presented in this paper have been run 
with two flavours of CR-injection, based on a high ef¬ 
ficiency model by Kang & Jones ( 20071) (model “0”) 
and on a lower efficiency model by Kang & Ryu|(|2013]) 
(model “1”). 

Radiative cooling (“c” in the name descriptor of 
Table [TJ has been modelled assuming a primordial 
composition of a fully ionized H-He plasma with a uni¬ 
form metallicity of Z = 0.3 Z 0 (where Z Q is the so¬ 
lar metallicity), using the APEC emission model (e.g. 
Smith et al. 2001). For the cold gas, with tempera- 
ture T V 1 0 4 K, we use the cooling curve of Smith 
et al. (|2011 [), which is derived from a complete set of 
of metals (up to an atomic number 30), obtained with 
the chem ical network of the photo-ionization software 
Cloudy (jFcrland et al. 1998 |). In order to model the UV 
re-ionization background (Haardt & Madau 1996), we 
enforced a temperature floor for the gas in the reash ift 
range 4 ^ ^ 7, as discussed in Vazza et al. (2010). 

All radiative simulations presented here aTso in¬ 
clude the effect of feedback from AGN, which 
are placed in high-density peaks (re ^ ttagn = 
10~ 2 cm~ 3 ). This method by-passes, both, the problem 
of monitoring the mass accretion rate onto the central 
black hole within each galaxy, and the complex small- 
scale physical processes which couple the energy from 


the black hole to the surrounding gas (i.e. the launch¬ 
ing of strong winds due to the radiation pressure of 
photons from the accretion disc). This is unavoidable, 
given that our best resolution is orders of magnitudes 
larger than the accretion disc region, let alone the dif¬ 
ficulty of modelling the radiative transfer of photons 
from the accretion region. At each feedback event we 
add an energy of £agn ~ 10 59 erg during a single 
time step, which typically raises the temperature in the 
cell to ~ 5 • 10 7 — 5 • 10 8 K at the injection burst. We 
model the feedback as a “bipolar thermal outflow”, i.e. 
the thermal energy is released into the ICM along two 
cells on opposite sides of the AGN cell, and the di¬ 
rection of the two jets is randomly selected along one 
of the three coordinate axes of the simulation. In the 
runs analysed for this work we include two slightly 
different implementations of AGN feedback in simula¬ 
tions with radiative cooling, which produce rather dif¬ 
ferent cluster scaling relations at z = 0. In the first run 
(“cl”), we trigger thermal feedback whenever the local 
maximum density exceeds the (comoving) gas density 
threshold, tiagn. starting from z = 1. In the second 
run, the AGN feedback is stalled at z ^ 4 (run “c2”). 
While the first approach is found to quench the radia¬ 
tive cooling inside most of halos but fails to produce 
cluster scaling relations that match observations (i.e. 
the M — T and Lx — T relations are much flatter than 
observations), the second produces a tilt in scaling re¬ 
lations that agrees well with observations, and is thus 
our preferred AGN-feedback method. 

At all redshifts the heating by the reionizing back¬ 
ground largely exceeds the heating contribution from 
CRs, which lose energy via interactions with thermal 
particles of the ICM (|Guo & Oh|2008j). This energy ex¬ 
change is modelled at run-time in our method for CRs, 
but it becomes important only for high gas densities 
( p/(/rmp) > 10~ 2 cm~ 3 ) tha t are typical of galactic 
cool cores ( Vazza et al.|2014b]). 


The runs analysed here assume a WMAP 7-year 
cosmology with Qq = 1.0, Qb = 0.0455, Qdm = 
0.2265, = 0.728, Hubble parameter h = 0.702, 

a normalisation for the primordial density power spec¬ 
trum eg = 0.81 and a spectral index of n s = 0.961 for 
the primordial spectrum of i nitial matter fluctuat ions, 
starting the runs at Zi n = 30 (Komatsu et a! .|201 1 i. 

Table [T] lists the simulations used for this work, 
showing the different choices for the physical mech¬ 
anisms, as well as the grid size. We simulated three in¬ 
dependent cosmological volumes with boxes of sides 
216 Mpc/h, 108 Mpc/h and 54 Mpc/h, respectively. 
In order to monitor resolution effects for each case, 
we produced re-simulations of each of these with a 
different number of cells (from 2048 3 to 256 3 ) and 
DM particles. In general, the missing effects of non- 
gravitational physics, such as radiative gas cooling and 
AGN feedback, are much more evident at high resolu¬ 
tions, where gas density and hence radiative losses are 
higher. For this reason, our study of the impact of non- 
gravitational physics is based largely on the smallest 
box (with side 54 Mpc/h). On the other hand, the im¬ 
pact of different CR injection efficiencies on the outer 
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cluster profiles is also captured at a lower spatial reso¬ 
lution, and therefore all our boxes include a study using 
different acceleration recipes. 


3 IDENTIFICATION OF FILAMENTS 


Our procedure to identify filaments uses the gas mass 
density to separate over- from under-dense regions. 
With this criterion, a single physical parameter deter¬ 
mines the identification procedure, with no further ar¬ 
bitrary settings of other physical quantities. The selec¬ 
tion is further refined by considering the numerical res¬ 
olution of the simulations and the expected geometry 
of a filament. These additional criteria mainly affect 
small objects (of typical size less than 1 Mpc 3 ). The re¬ 
sulting methodology is implemented based on a combi¬ 
nation of different algorithms developed and optimised 
in the context of material interface reconstruction. The 
accurate identification of the interface between two or 
more materials represents a major issue in these appli¬ 
cations since a single computational cell can be com¬ 
posed of pieces of different materials. Various solu¬ 
tions have been developed in order to reconstruct in¬ 
terfaces, both for simulations and for visualisation and 
data analysis. Examples are pr ovided by Meredith & 
Childs] <[2010]), |Fujishiro et al.| <J1995]>, |Hamson et al. 
(|2011 1 , Anderson et al.|(]2(TTU I and references therein. 


We can identify two distinct material phases: the 
warm/hot gas which fills collapsed cosmological struc¬ 
tures on one side, and the cold under-dense phase typ¬ 
ical of voids on the other side. The accurate represen¬ 
tation of the boundaries between these two phases is 
the key to define the properties of the identified ob¬ 
jects. In this work, we have adopted an Isovolume- 
based approach, described in |Meredith (2004). A quan¬ 
titative evaluation of this approach and of its perfor¬ 
mance, compared to other methodologies, can be found 
in|Meredith & Childs|(j2010]). It has been further eval¬ 
uated for our specific needs. The resulti ng a ccuracy 
and performance are presented in Section 13)21 and Ap¬ 
pendix A respectively. 

In order to implement a computationally efficient, 
flexible and extensible filament reconstruction proce¬ 
dure, we have exploited the Visit data visualisation 


and analysis framework (Childs et al. 2011). Visit is 
designed for the efficient processing of large datasets 
through a combination of optimized algorithms, the 
exploitation of high-performance computing architec¬ 
tures, in particular through parallel processing, based 
on the MPI standard, and the support of client-server 
capabilities, allowing efficient remote visualization. It 
can be used interactively, through a graphical user in¬ 
terface, and it can be scripted using the Python pro¬ 
gramming language in order to automate the data pro¬ 
cessing. Finally, it natively supports a number of file 
format adopted by popular astrophysical simulation 
codes, such as Cosmos, FLASH, Gadget, Chombo (be¬ 
sides, of course, ENZO ) and in general all files adopt¬ 


ing the VTK format. [j]We refer to the official documen 
tation[£]for all the technical details on the software. 

Visit includes the aforementioned Isovolume algo 
rithm, which, combined with the Connected Compo¬ 


nents _filter (also provided by the software, Harrison 


et al.||20Tij), allows the segmentation of the data into 


distinct objects. An example of a filament identified 
and reconstructed by our procedure is shown by Fig¬ 
ure [Q 


3.1 Filament identification procedure 

The main steps performed by our implemented identi¬ 
fication procedure can be summarised as follows: 


(i) Identification of large clumps. In a first step we 
separate filaments from galaxy clusters. Clusters can be 
described as high-density clumps in the matter distri¬ 
bution. Such clumps can be identified first by selecting 
isovolumes with 


£BM . 

ft ®cll 

00 


(1) 


where pgM is the cell’s baryon mass density, a c i is a 
proper threshold and po is the critical density at present 
time. For each clump, Visit returns the volume V c \ 
and the coordinates of its centre. At this point, assum¬ 
ing spherical symmetry, we calculate the radius of each 
cluster as: 


R = 



1/3 


( 2 ) 


where /5 is a multiplicative factor needed to rescale the 
radius to values of the order of 1 Mpc (expected for 
galaxy clusters). Setting a c i = 100 and ft = 10 proved 
to be effective for a proper cluster characterization. 

Clusters with radius R ft 1 Mpc are discarded by 
removing all cells falling within spheres with centre in 
the cluster centroid and radius R. 

(ii) Filament identification. The Isovolume algo¬ 
rithm is used once more on the residual cells to iden¬ 
tify the volumes that obey Equation [T] only with a 
new threshold am- The exact results of the identifi¬ 
cation procedure depends critically on thi s pa rameter, 
whose setting will be discussed in Section|33]and Ap¬ 
pendix A. The connected components algorithm is then 
used to combine cells belonging to distinct (i.e. non in¬ 
tersecting) filaments, assigning to each filament an Id 
(an integer number) and marking each cell belonging 
to a filament with that Id. For each filament, the vol¬ 
ume (Vfii), the mass (Mm) and the average temperature 
(7m), are calculated. 

(iii) Small clumps removal. The next step consists 
of removing small objects with volumes Vm < V res . 
In order to have a unique volume cutoff for all sim¬ 
ulations, the parameter V r res is set equal to the lowest 
spatial resolution of our runs, V/ es =(0.3 Mpc) 3 . This 


1 http://www.visitusers.org/index.php?title=Detailed_list_of_file_formats_VisIt_supports 

2 https://wci.llnl.gov/codes/visit/manuals.html 
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Filaments in Cosmological Simulations 5 


Table 1 . List of the simulations run for this project. Column 1: box size. Column 2: number of grid cells. Column 3: spatial resolution. Column 4: physical 
implementation. Column 5: identification name of each run. The ID must be read as follows: the first number refers to the box size, 1 referring to 216 Mpc/h, 
2 to 108 Mpc/h, 3 to 54 Mpc/h. The second symbol refers to the physical model. Finally, the third number indicates the ID computational mesh size in cells. 
These IDs will be adopted throughout the paper. 


ibox [Mpc/h] 

A^grid 

Ax[kpc/h] 

physics 

ID 

216 

2048 3 

105 

non-rad.+CR(KR 13) 

1 -1.2048 

216 

1024 3 

210 

non-rad.+CR(KR 13) 

1-1.1024 

216 

1024 3 

210 

non-rad.+CR(KJ07 ) 

1-0.1024 

108 

1024 3 

105 

non-rad.+CR(KR13) 

2-1.1024 

108 

1024 3 

105 

non-rad.+CR(KJ07) 

2-0.1024 

108 

1024 3 

105 

cool 1 ,+AGN feedback+CR(KR 13) 

2-c 1.1024 

108 

1024 3 

105 

cool2.+AGN feedback+CR(KR 13) 

2-c2_1024 

54 

1024 3 

52 

non-rad.+CR(KR 13) 

3-1.1024 

54 

1024 3 

52 

non-rad.+CR(KJ07) 

3-0.1024 

54 

1024 3 

52 

cool 1 ,+AGN feedback+CR(KR 13) 

3-cl_1024 

54 

512 3 

105 

non-rad.+CR(KR 13) 

3-1.512 

54 

512 3 

105 

non-rad.+CR(KJ07 ) 

3-0.512 

54 

512 3 

105 

cooll.+AGN feedback+CR(KR13)+AGN 

3-c2_512 

54 

256 3 

210 

non-rad.+CR(KR 13) 

3-1.256 



Figure 1 . Example of a filament identified by the Isovolume method. Part of the object is outlined by the reconstructed surface mesh. The rest of the object is 
coloured with the temperature. 
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ensures that a filament is always larger than a single 
cell. 

(iv) Shape selection. Although the previous steps 
arc designed to eliminate small clumps and large 
galaxy clusters, nonetheless a few outsiders can still 
be present at this stage. Hence, we use the following 
two-stage cleaning procedure to remove all remaining 
round-shaped, isolated structures. First, we calculate 
the bounding box enclosing each identified object and 
accept as filaments all those objects with: 

MAX(r X j, r IZ , r^) > cr, (3) 

where r a b is the ratio between axes a and b of the 
bounding box. A value of a = 2 ensures that the se¬ 
lected object stretches in a specific direction without 
being too restrictive. The second stage is intended to 
recover all those objects which do not match criterion 
purely based on axis ratio. This happens, for instance, 
to filaments laying along the bounding box diagonal. 
The second shape selection consists of calculating the 
filling factor Fy of each remaining object, defined as: 

Fy = (4) 

where Vb ox is the bounding box volume. All objects 
with Fy < ip arc classified as filaments. The value of 
c p is set according to the ratio between the volume of 
a cubic bounding box and that of a cylindrical filament 
lying on its diagonal, with a base radius set as a frac¬ 
tion, namely 1/5, of the box side. The values of the 
parameters a and p are discussed in Appendix A. 

(v) Data export. The properties of the cells belong¬ 
ing to each filament, such as the mass density, tem¬ 
perature, velocity, geometric and topological informa¬ 
tion, are exported to output files using the VTK for- 
ma{3 Since the number of cells identified as filaments 
is some fraction of the total computational mesh, the 
file size is reduced correspondingly, making further 
data management and processing significantly easier 
and faster, which is very convenient for the large runs 
we use in this work. 


3.2 Validation 

We validate the Isovolume algorithm using spherically 
symmetric mass distributions modelled by Gaussian 
functions, representing idealised clusters and clumps, 
and cylindrical shapes to mimic filaments. For the fil¬ 
aments, the mass distribution orthogonal to the main 
axis is again modelled as Gaussian. We have first eval¬ 
uated the accuracy of surface areas, volumes and den¬ 
sity estimates as a function of the mesh resolution for 
single objects. Table [2] shows the error, calculated as 
the fractional difference between the calculated and the 
exact value of each of the three quantities, as a func¬ 
tion of the radius, expressed in computational cells. 
For both spheres and cylinders, the density estimate is 
always accurate to ^ 1%. Volumes and surface areas 


3 http://www.vtk.org/VTK/img/file-formats.pdf 


Table 2. Fractional error in the calculation of surface areas (columns 2 and 
5), volumes (column 3 and 6) and densities (column 4 and 7) for spheres 
and cylinders as a function of the distance (in cells) from the centre/axis 
(column 1) of the sphere/cylinder. 


R 

<5 A/A 

sv/v 

Sq/q 

5A/A 

SV/V 

Sq/q 


Sphere 



Cylinder 



3.84 

0.0762 

0.1203 

0.0024 

0.0616 

0.0836 

0.0010 

6.40 

0.0265 

0.0421 

0.0011 

0.0423 

0.0453 

0.0007 

12.80 

0.0049 

0.0079 

0.0002 

0.0306 

0.0183 

0.0010 

25.60 

0.0003 

0.0004 

0.0002 

0.0200 

0.0048 

0.0069 

38.40 

0.0014 

0.0020 

0.0021 

0.0067 

0.0207 

0.0193 


are more sensitive to the resolution (see also Meredith 
|& Childs||201 Oj ). For radii smaller than 4 cells, the er- 
ror approaches (or is even larger, as in the case of the 
sphere’s area) 10%. 

Datasets with various combinations of spheres and 
cylinders have been used to model the matter distri¬ 
bution from a cosmological simulation. For all tests, a 


computational mesh of 128 3 cells is used. An example 
is presented in Figure [2j where isovolumes at agi = 1 
and dfii = 50 are shown in the left and right panels, 
respectively. The accuracies of the volume and density 
estimates have been calculated as a function of agi, ob¬ 
taining the results presented in Table [3] Volume and 
density estimates are compared to those obtained by a 
simple algorithm which sums the contributions of all 
the cells above the given threshold, run at much higher 


resolution (1024 3 ). In all tests, the number of identified 
objects at each density threshold is correct. Only ob¬ 
jects crossing the faces of the periodic computational 
box could not be counted properly because our pro¬ 
cedure does not yet support periodic boundary condi¬ 
tions. In this case, a single object is split into multiple 
components, each treated as a separate filament. This 
leads also to underestimate the length of any filament 
crossing one of the box’s side, as well as to the artifi¬ 
cial increase of small-scale objects because of this ar¬ 
tificial fragmentation. However, this problem is statis¬ 
tically not very relevant, as our volumes are typically 
much larger than the largest objects in the volume. 

The differences in the volume estimates are be¬ 
tween the 1.6 and 4.7%. The maximum difference is 
at the highest density thresholds where the size of the 
various objects approaches the mesh resolution. The 
minimum is for ag] = 5.0 with objects that are still 
extended and have a simple geometry (disjoint spheres 
or cylinders). Density estimates show differences of the 
order of or less than 2%. 


3.3 Density threshold and spatial resolution 

The parameter ag] is critical for the characterization 
of the filament distribution. Its value is discussed be- 
low with additional considerations presented in Section 
4.3|and in Appendix A. This parameter determines, to¬ 
gether with the numerical spatial distribution, the fea¬ 
tures and the connectivity of the detected structures. 
Figure 3] shows the mass distribution obtained using 
four different values of ag], ranging from 0.5 to 2 for 
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Figure 2. Isovolumes at a^\ = 1 (left) and a^\ = 50 (right) extracted from a test dataset combining spheres and cylinders with different mass distribution 


Table 3. Fractional error in the calculation of volumes (column 2) and den¬ 
sities (column 3) for combination of spheres and cylinders as a function of 
afii (column 1). 


“fil 

sv/v 

Sq/q 

0.25 

0.026 

0.010 

0.50 

0.025 

0.011 

0.75 

0.024 

0.012 

1.0 

0.023 

0.013 

2.0 

0.020 

0.018 

5.0 

0.016 

0.022 

50.0 

0.047 

0.004 


the model 2-1.1024, keeping all else constant. Similar 
results are obtained for the other models (not shown). 
Values of the mass density wit hin this range are those 


expected for filaments (see e.g. Cen & Ostriker 1999 
Dave et al.|20011 >. 


A qualitative visual inspection of the distributions 
suggests the adoption of a value of am in the range 
0.75-1.0, smaller values leading to structures perco¬ 
lated across the whole computational box, and higher 
values generating clumpy distributions with limited (or 
even no) connectivity. 

Figure [4] shows the number of objects extracted 
from the simulations with box size of (50Mpc/h) 3 as 
a function of am- Results obtained with computational 


meshes of 1024 3 , 512 3 and 256 3 cells are presented, 
in order to investigate also the impact of spatial reso¬ 
lution on the results. Spatial resolution directly affects 
the number of objects identified by our algorithm, the 
differences decreasing at large values of the threshold 
parameter afu. The number of filaments increases by a 


factor of ~ 3.5 for am = 0.5 going from the 256 3 to the 
1024 3 run, while the increase is only by a factor ~ 2 
for the largest threshold, am = 5. This shows that reso¬ 
lution strongly affects the identification of low-density 
filamentary structures below a given spatial scale, iden¬ 
tifiable at about 0.2 Mpc, due to the larger diffusion of 
the gas at that resolution, which smoothens the mass 
density in the outer parts of the filaments, connecting 
structures otherwise distinct and leading to an overall 
decrease in the number of detected objects. As a con¬ 
sequence, among our datasets (Table 1), we expect the 
resolution to affect mainly the coarsest runs, with spa¬ 
tial resolution of 210 kpc/h (1-0.1024, 1-1.1024 and 
3-1.256). 

It is worth noting that for values of am between 
0.5 and 1, the various physical processes acting in the 
different simulations have the smallest influence on the 
statistics. Furthermore, at the highest resolution (1024 3 
cells) the value am = 1 represents an inflection point 
for all the runs, suggesting a change in the properties 
of the distributions under investigation. 

Figure [5] (left panel), shows the dependence of the 
average mass density of a filament, p avg , on the pa¬ 
rameter am- As expected, the quantity o avg grows with 
increasing am since lower-density regions are progres¬ 
sively cut out. The growth is almost linear and the slope 
higher for the “cl” model, due to the cooling processes 
which dominate the AGN feedback. In model “c2”, 
that has effective AGN feedback, the effect of cool¬ 
ing is completely compensated by the energy injection, 
which raises the pressure and stops the infall of matter. 
The runs with CRs, “0” and “1”, give similar results 
at all resolutions. This indicates that CRs only have 
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Figure 3. Mass distribution obtained using clr\ = 0.5 (top left panel), 0.75 (top right), 1.0 (bottom left) and 2.0 (bottom right), for the model 2-1.1024. Each 
different object identified by our algorithm has a different colour. 


a small effect on the dynamics of the gas within fila¬ 
ments. The model “cl” has densities higher than the 
other models, due to the much lower thermal support 
within the collapsing structures, leading to highly com¬ 
pressed matter distributions. The trend with resolution 
of the mean gas density shows that for most of the in¬ 
vestigated threshold values am (with the exception of 
am > 4 extractions, which should be dominated by 


low-number statistics) the mean density increases with 
resolution, due to the higher achievable compression 
and the increased presence of substructures within fil¬ 
aments. 

The right panel of Figure|5]shows the average tem¬ 
perature of gas in filaments, Ta Vg , as a function of am- 
Also the temperature increases with the density thresh¬ 
old parameter, faster for am < 1.0. At higher val- 
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a fii 


Figure 4. Number of filaments per 100 cubic Mpc as a function of a^\ 



a fil 



Figure 5. Average baryonic mass density (left panel) and temperature (right panel) in filaments as a function of a^\. 


ues, the temperature grows slowly with density. Inner 
parts of the filaments are thermalized at an almost con¬ 
stant temperature by shock waves propagating from the 
centre outwards. The outer part of the filaments, cap¬ 
tured only at low afii, can comprise non-shocked cells 
that lower the average temperature. The most striking 
feature of the graph is the trend of temperature with 


resolution. The lowest resolution run (256 3 ) presents 
the largest average temperature, at all investigated am, 
and further increase in resolution is followed by a de¬ 
creased mean temperature. This trend is caused by the 
weakening of outer accretion shocks (yielding a lower 
thermalisation efficiency in filaments) when resolution 


is increased, an effect already mentioned in Vazza et al. 
(|2014b|). On the opposite end, the model 3-c IJ024 has 


the lowest average temperatures at all values of p avg , 
as a result of the strong effect of cooling. The “c2” 
model has a peculiar behaviour. At low densities, cool¬ 
ing tends to dominate AGN feedback and the resulting 
average temperatures are clearly below all the other 
models (all but “cl”). However, when higher-density 
objects are selected, AGN heating dominates, leading 
to the highest temperatures (neglecting the 256 3 case). 
The analysis of these results indicates that overall a 
value afii = 1 is a suitable choice to identify filamen¬ 
tary structures. This choice of parameter is further sup¬ 
ported by the analysi s of the density and mass profiles 
discussed in Section |4.3| Hence, this will be the fidu¬ 
cial value for the rest of the paper. 
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Figure 6. Eight filaments selected from simulation 3-1-1024. Each object is identified by a different ID (an integer number). The top left panel shows 
the selected filaments altogether in the simulation box, for comparison. The other panels are close-ups on each single filament. Colours represent the gas 
temperature in the outer surface of each object. 


4 PROPERTIES OF THE FILAMENTS 

By keeping our fiducial density threshold of agi = 1.0 
we identified all filaments in the simulated datasets and 
computed their average properties, as a function of res¬ 
olution and of the adopted physics. About 30000 fila¬ 
ments were identified in our largest run (1-132048) and 
about 1000 objects in our highest resolution runs. 


4.1 Visual analysis 

In Figure [6} we show eight objects extracted from one 
of our highest resolution runs (3-1 _1024), as a repre¬ 
sentative sample of the filaments identified by our pro¬ 
cedure. The eight objects are between 10 and 28 Mpc 
long, the longest being that with ID=7887 (27.9 Mpc). 
In the top left panel the selected filaments are presented 


together in order to show their location and to compare 
them to each other. The other panels are close-ups on 
each single object. Colours represent the temperature 
of the density isovolumes. The images show how fila¬ 
ments can have heterogeneous geometries, depending 
on the environment in which they lie. In some cases, as 
for filaments 678, 6397, and 2957, they have elongated 
and rather regular shapes, with comparable thickness 
and surface temperature. They have one or two main 
blob-like structures and they typically bridge pairs of 
massive galaxy clusters. We indicate such structures 
as “giant bridges”. On the other hand, filaments 7887, 
3153 and 14599 have a more branchy structure, with 
multiple blob-like components and very thin segments, 
which present a systematically colder surface, and sev¬ 
eral branches. Another peculiar type of filament is rep¬ 
resented by filament 2131, which is more than 16 Mpc 
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Mass 


Figure 7. Relation between the enclosed gas mass and volume for the filaments in our sample (only 1/20 of filaments are shown for each model, for a clearer 
view). The additional lines in colour show the best-fit relation within each sample, while the two grey lines show the best-fit for the population of galaxy 
clusters extracted from run 2-l_1024 and 2-cl_1024. 




Figure 8. Relation between the enclosed gas mass and the estimated length of filaments in our sample. The left panel show the full range of values in the 
dataset (only 1/20 of filaments are shown for each model, for a clearer view). The additional lines in colour show the best-fit relation within each sample. The 
right panel only focuses on the objects with mass M ^ 1O 13 M0 (in this case all objects are plotted, same colour coding as in the left panel). The length of 
each object is an estimate based on the bounding box of each filament, as described in Section[3] 


long and is much thinner and colder than all the other 
objects. The nature of these objects require specific 
investigation to be fully understood. However, Their 
features point towards filaments in a early evolution¬ 
ary stage, with ongoing mass accretion, heated up by a 


first generation of weak shocks propagating in a rather 
under-dense environment, and a total enclosed over¬ 
density which is just above our detection threshold. 
This class of filaments can be called “long-thin”. Fi¬ 
nally, as in the case of filament 8975, we find objects 
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Figure 9. Relation between the enclosed gas mass and the estimated redius of filaments in our sample. The left panel show the full range of values in the 
dataset (only 1/20 of filaments are shown for each model, for a clearer view). The additional lines in colour show the best-fit relation within each sample. The 
right panel only focuses on the objects with mass M ^ 10 13 M® (in this case all objects are plotted, same colour coding as in the left panel). The fit for 
the 3-c 1.1024 model is omitted since statistically not meaningful (too few points available in the selected mass range). The radius of each object is estimated 
assuming cylindrical symmetry, as Rfu = (V/irL) 1 ' 2 , where V and L are the volume and the length of the filament respectively. 



Mass 


Figure 10. Relation between the enclosed gas mass and the average temperature for the filaments in our sample (only 1/20 of filaments are shown for each 
model, for a clearer view). The additional lines in colour show the best-fit relation within each sample, while the two grey lines show the best-fit for the 
population of galaxy clusters extracted from run 2-1.1024 and 2-c 1.1024. 
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for which it is difficult to recognise the typical charac¬ 
teristics expected for filaments. They are found in rich 
regions and connect multiple larger halos. As for the 
other classes, they contain large almost-spherical mat¬ 
ter clumps. Objects with these features are classified as 
“irregular”. 

Irrespective of the classification, most of our fil¬ 
aments contains spherical clumps that have low av¬ 
erage over-densities (around 20-50) and temperatures 
(around 0.1 keV). Furthermore, they have radii smaller 


than 1 Mpc and masses around 10 12 M 0 or lower. They 
can thereby be considered as proper sub-structures of 
filaments, at least at the resolution we can achieve. 

According to the above simple classification, we 
have estimated that about 30% of the filaments in our 
simulation box can be classified as “giant bridges”, 
while approximately the 36% are “branchy”. “Long- 
thin” objects are less frequent, but still account for the 
12% of the filaments. “Irregular” objects represent the 
16% of the population. It must be stressed that this 
classification is performed visually on a limited sam¬ 
ple of large objects (about 100 objects longer than 7 
Mpc), while smaller objects tend to be more difficult 
to classify, as resolution effects coarsen their substruc¬ 
tures and smooth out their shape. 

The remaining objects comprise isolated clumps 
emerging from the gravitational collapse process and 
spurious leftovers of the filament reconstruction, in 
part icular at the edges of our boxes (see also Section 


4~2| ). This also affects the low-mass end of the M-T re¬ 
lation as will be discussed in the following sections. 


4.2 Statistical properties of the filaments 

Figure [7] shows the distribution of the gas mass of fila¬ 
ments versus volume. Most objects fall along a narrow 
scaling relation, fitted by V ~ 1000 Mpc 3 ■ M/(2 ■ 
10 14 Af©). For comparison, we also show the best fit 
relations for the case of galaxy clusters, extracted from 
the 2-1.1024 and 2-c 1.1024 runs (for a total of ~ 300 
objects in the gas-mass range 10 12 ^ M/Mq ^ 10 15 ). 
Independent of resolution and physics, a self-similar 
scaling relation for filaments is found, as for galaxy 
clusters, yet with a higher normalisation. However, the 
normalisation of the best fit is smaller for the runs with 
efficient cooling (“cl” models), implying that cool¬ 
ing produces significantly more compact filaments, and 
that this removes more volume from the WHIM phase 
and locks it into overcooling halos, in turn reducing the 
volume of the filaments. 

Figure [8] gives the distribution of gas mass versus 
the length of each filament, which is approximated us¬ 
ing the diagonal of the bounding box defined by the 
reconstruction algorithm (Section [3]). The distribution 
for the full sample of each run is also described by 
a well-defined power-law with a rather narrow disper¬ 
sion of values around the mean, and with a self-similar 
scaling with mass, L oc M 1 / 3 . The presence of unre¬ 
solved small-size halos embedded within irregular fil¬ 
aments can artificially steepen the relation by increas¬ 
ing the total mass of gas in filaments, while for higher 


masses their contribution becomes negligible. There¬ 
fore, we also show the distribution and best fit (which 
become just a bit steeper, but still consistent with the 
oc M 1 / 3 scaling) for the distribution of filaments lim¬ 
ited to ^ 10 13 AT., objects, which are characterised by 
L ^ 10 Mpc. However, notice that our method system¬ 
atically underestimates lengths, hence this trend repre¬ 
sent a lower bound for the mass-length relation. The 
runs with cooling show a rather steeper scaling, sug¬ 
gesting that filaments of the same length are less mas¬ 
sive than in non-radiative runs. This is easily explained 
by the mass drop-out of cooling gas onto clumps within 
filaments, which can significantly remove a fraction 
of the WHIM from the most diffuse phase, and at the 
same time promote the full collapse of small-size halos 
and their detection by our algorithm. 

The radii of the filaments, presented in Figure [9j 
arc estimated assuming cylindrical symmetry, such that 
R = {V/-kL ) 1 / 2 , where R, V and L are the radius, 
the volume and the length of the filament, respectively. 
This is a rough estimate since many objects are far 
from a simple cylindrical geometry. However, it pro¬ 
vides an indication of the characteristic transversal size 
of a filament and its distribution with mass, with radii 
between 0.1 and 0.3 Mpc for objects up to 10 12 AT.., 
growing to more than 1 Mpc for the largest filaments, 
with masses between 10 14 and 1O 15 M 0 . If resolution is 
increased, the radii of filaments shrink, as indicated by 
comparing the 1-0.1024 and 1-1.1024 to the 1-1.2048 
interpolation fits. This affects in particular' smaller ob¬ 
jects, as shown by the right panel of Figure [9| where 
we present the radius-mass relation for the sub-set of 
filaments with masses above 1O 13 M0. The transverse 
size of these large filaments is less affected by reso¬ 
lution. This indicates that the properties of those ob¬ 
jects whose radius is comparable to the cell size (low 
mass filaments) are not fully converged with spatial 
resolution. Also radiative cooling (run cl) leads to nar¬ 
rower filaments, with the smallest radii associated with 
3-cl_1024 and l-cl_1024. This effect is compensated 
by the efficient AGN feedback in the c2 models. 

For the same objects, we show in Figure [TO] the re¬ 
lation between the enclosed gas mass and the volume- 
averaged gas temperature. Here a less defined scaling 
relation is seen. The overall distribution follows a scal¬ 
ing law similar to galaxy clusters in the non-radiative 
case (T ~ M 2 / 3 ), with 3. rN_/ 1 order of magnitude 
lower normalisation. However, a small fraction of the 
identified objects has unexpected properties, with low 
masses (below 10 12 AT.,) and high temperatures above 
10 5 K). These outliers are due to our reconstruction pro¬ 
cedure, which does not support periodic boundary con¬ 
ditions. They correspond to small leftovers of larger 
(and hotter) objects crossing one of the sides of the 
computational box and there cut in two distinct parts. 

Variations of CR physics produce extremely small 
effects, and the scaling relations are essentially the 
same for the “1” and “0” models. Instead, the role 
played by gas cooling and feedback from AGN be¬ 
comes much more relevant compared to the previous 
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Figure 11. Number of filaments normalised to 100 cubic Mpc as a function of their volume. 



10 11 12 13 14 15 16 

log 10 M[M o ] 


Figure 12. Number of filaments normalised to 100 cubic Mpc as a function of their gas mass 
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Figure 13. Number of filaments normalised to 100 cubic Mpc as a function of average temperature 
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Figure 14. Number of filaments normalised to 100 cubic Mpc as a function of the pressure ratio between CR energy and gas energy. 
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scaling relations since temperature is directly affected. 
As for galaxy clusters, runs including cooling and 
AGN feedback produce a steepening of the scaling 
relation in the low-mass end of the distribution. Due 
to cooling, the average temperature at all masses is 
reduced significantly, and also the best-fit relations 
for the runs including cooling are steeper than in the 
non-radiative case. The effect is even bigger than in 
the case of galaxy clusters, as shown by Figure |T0} 
comparing the best fit for clusters and filaments in 
the same volume. This can be explained considering 
that first, below T ^ 10 5 K the radiative cooling is 
strongly affected also by line cooling and second that 
although AGN feedback is effective within clusters 
and can reach out to the WF1IM in filaments, it is ob¬ 
viously less efficient there given the distance (several 
~ Mpc) from its release. As also pointed out in the 
next Section, the cooling efficiency in this regime is 
probably overestimated in our case, since we assumed 
a (constant) metallicity of 0.3 Z@ everywhere in the 
simulated volume. Furthermore, the mass dropout into 
forming stars is not included. 


The volume/mass/temperature number distribu¬ 
tions for all objects are given in Figures 11 13 Apart 
from a deficit of objects with V ^ 1 Mpcr in the 
1-0.1024 and 1-1,1 024 runs (due to coarse resolu¬ 
tion effects, Section 13.31) , the volume distribution of 
objects is very similar m all runs and is described 
by a simple logio(Ngi) oc V -1 relation across ~ 3 
decades in volume (Figure m The various curves 
show a good numerical convergency down to gas 


masses of ~ 10 11 Mq indicating that our reconstructed 
distribution of objects is complete down to this mass. 
The total volume occupied by filaments ranges from 
3.8 to 4.5% of the computational volume in all our 
runs. 


When physical resimulations of the same vol¬ 
ume are compared, we we find a significant excess 
of large-volume filaments when the more efficient 
“c2” feedback model is used, for ^ 100 Mpc 3 as 
well as a significant deficit of halos in the range 10 
Mpc 3 ^ V ^ 100 Mpc 3 . This is due to the powerful 
AGN feedback at high redshifts which leads to the ex¬ 
pulsion of gas from the centre of dense proto-clusters 
and affects the properties of a few tens of massive 
filaments around the most massive halos in the box. 
The efficient feedback releases hot gas beyond the 
radius achieved by non-radiative runs, with a strong 
impact on the densest regions, where the AGN feed¬ 
back mostly takes place. This enriches even rarefied 
environments with additional gas, and leads to the 
junction of otherwise disconnected structures through 
the extra gas. A more detailed analysis will be given in 
the next section (e.g. Figure p~6j). 

On the other hand, the effect of radiative cooling is 
enough to produce more compact objects within fila¬ 
ments than in the non-radiative case, thereby causing 
an increase of the low-volume objects in the sample. 
Both effects are relevant for any filament detection 


scheme based on gas. Very similar trends are found 
for the distribution of filament mass (Figure [T2]). 
We conclude that on average only one filament with 
M ^ 1O 14 M 0 and V ^ 1000 Mpc 3 is found within a 
volume of (100 Mpc) 3 , while within the same volume 
we can find ~ 10 objects with M ^ 1O 13 M 0 and 
~ 10 2 objects with M ^ 1O 12 M 0 . 


The distribution of temperatures (Figure 131 
clearly shows the imprint of radiative cooling on me 
WF1IM. Cooling removes a significant part of the cos¬ 
mic population of filaments from the soft X-ray bands 
and pushes it much below 10 5 — 10 6 K (still above 
the minimum level set by our assumed UV heating 
background, which had dropped off at a level of a few 
~ 10 2 K at z = 0), as already seen in the average tem¬ 
peratures in Figure [5] The coarse resolution in the 1- 
0.1024 and l-l_102Truns is responsible for an excess 
of large temperatures compared to the other runs with 
a higher resolution. 

Finally, we calculate the ratio between the thermal gas 
pressure and the pressure of CR-protons accelerated by 
shock waves, as in Vazza et al. (j2014b|. FigurefT4|gives 
the number distribution of filaments and compares all 
our “1” and “0” models for the non-radiative case. For 
the bulk of the population the pressure ratio within the 
volume of filament is Pqr/P g ~ 0.1 — 0.2, with a ten¬ 
dency of larger filaments to have lower ratios. As ex¬ 
pected, the lower efficiency “1” model (|Kang & Ryu 


20131 produces a smaller CR-pressure compared to 
the higher efficiency “0” model (|Kang & Jones|2007j), 
Flowever, the differences are small. At the same reso¬ 
lution, the two distributions differ by ~ 20 — 30 per¬ 
cent. This is expected because in the different simula¬ 
tions the enrichment of CRs into filaments is equally 
strong, owing to the fact that the saturated efficiency of 
the two acceleration models is only different by a fac¬ 
tor ~ 30 percent in our treatment (Vazza et al.j2014bj). 
The effect of increasing resolution is to lower the av¬ 
erage pressure ratio due to the overall weakening of 
shocks, which is made manifest by the shift of the peak 
of this relation to lower values, going from the coarse 
1-0,1024 and 1-1,1024 runs (where the number of ob¬ 
jects is also found to be smaller than in the other runs) 
to the 3-0,1024, 3-1,1024 and 3-c 1,1024 runs. This 
suggests that very large values of the pressure ratio are 
mostly driven by resolution effects on the shock pop¬ 
ulation 

thorough investigations 
tribution of CRs in these runs have been already pre¬ 
sented Vazza et al. (2014b[, Figure 11) and showed that 
Pcr/P g 1S a rather well converged quantity (within a 
factor ~ 2 — 3) for all across the full range of cosmic 
environment, for spatial resolution equal or better than 
~ 200 —300kpc, i.e. well in the range of what explored 
by the runs we use in this work. Overall, these findings 
are in agreement with Vazza et al. (20l4b), and indi¬ 
cate that, provided diffusive shock acceleration (DSA) 
can take place at low overdensities, cosmic filaments 
can store a large amount of CRs. 


(|Ryu et al.|2003[ Vazza et al .|2011j >. Flowever, 
hi investigations on resolution effects in the dis- 
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Figure 15. Phase diagram ( p,T ) for the filaments FI (top panels) and F2 (bottom panels), comparing the non-radiative run and in the cooling and feedback 
run (c2). The colour coding gives the fraction of volume characterised by a given phase. For each object we also show the contours of the phase diagram for a 
larger « (35 Mpc) 3 box in order to compare with the entire distribution of cosmic baryons. 


4.3 Properties of individual filaments 


In order to monitor the impact of non-gravitational 
physics (cooling and feedback from AGN) onto fila¬ 
ments, we analysed the properties of individual objects. 
For this purpose, we extracted from the 3-1 1024 and 
3-cl_1024 runs all cells belonging to two objects: fila¬ 
ment “FI” with a length of ~ 15.5 Mpc and filament 
“F2” ~ 10.2 Mpc long. The maps of projected gas tem¬ 
peratures for these objects are shown in Figure [161 For 
the two different physical models, the filamentsrnor- 
phologies are quite similar. However, the combined ef¬ 
fect of radiative cooling and feedback produces more 
substructures within each filament which are narrower 
compared to the non-radiative runs. Also the presence 
of more compact sub-halos causes a higher fragmenta¬ 
tion in several portions of the filaments in the radiative 
runs. 


Other thermodynamical differences are more con¬ 
spicuous in phase diag rams, as in the ( T,p ) phase di¬ 
agrams of Figure [15] The properties of the cells in 
filaments can be compared to the global phase dia¬ 
gram of the larger « (35 Mpc) 3 boxes containing the 
two filaments, in order to compare with the phase di¬ 
agram of cosmic baryons in the volume. The two ob¬ 
jects show a very similar thermodynamic pattern and 
contain most of the dense and warm-hot intergalactic 
medium (T ^ 10' K) outside galaxy clusters. Most 
of the volume in these two objects is taken up by gas 
with densities p ~ 10(p) and temperature T 10 6 K. 
The combined action of cooling and feedback is rele¬ 
vant in both cases. Cooling promotes the formation of 
denser clumps in both filaments (visible as horizontal 
stripes on the phase diagrams), and it is found to af¬ 
fect the filaments even at T C 10' K, i.e. in a regime 
where metal cooling can be dominant and that future 
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Figure 16. Maps of projected temperature for filament FI and F2, for both 
non-radiative (left) and cooling & feedback (right) runs. The white con¬ 
tours give the boundaries of the filaments as reconstructed by our algo¬ 
rithm, while the white points show the spine of each filament, that we used 
to compute the profiles in Figure p~7|19| Each panel has axes is cell units 
(here Ax = 74 kpc). 


X-ray telescopes (e.g. ATHENA, Nandra et al. 2013 > 
will probe. Although recent observations ot metaiiicity 
at the edge of the Perseus cluster found a mean metai¬ 
iicity of the order of ~ 0.2—0.3 Zq ([Urban et al.|2014|), 
our assumption is likely to overestimate the presence of 
metals at early times, and translates into an upper limit 
on the amount of radiative cooling of the WHIM (e.g. 


Cen & Ostriker||2006[ ISmith et al.||2011t |Maio et al. 
201 ip . 


At the same time, the effect of feedback from the 
nearby clusters is enough to increase the temperature 
by a factor of 2-3 in most of the volume compared 
to the non-radiative run. The most peripheral parts of 
some clusters expand due to the AGN outbursts and, 
by getting more rarefied, are identified as part of fila¬ 
ments by our algorithm (this is particularly evident in 
the upper part of filament F2). This leads to a slight in¬ 
crease of the size of the structure due to percolation of 
expanding shells of matter, which in turn enhances the 
maximum length of several filame nts i n our radiative 
runs (as also pointed out in Section[4~2]>. 


FI & F2 



Figure 18. Transverse profile of the enclosed gas mass for the filaments 
FI (solid) and F2 (dashed), in the non-radiative run and in the cooling and 
feedback run (c2). The grey arrows show the approximate location where 
the mean density profile of filaments is equal to the different values of the 
density thresholds a^i- The profile of each filament has been normalised for 
the total enclosed mass to the last radial bin. 


The effect of our prescription for the thermal feed¬ 
back and the fact it reaches out to the outer part of fil¬ 
aments is similar to what is reported by other authors,_ 
who investigated similar fee dback schemes (Kang et al. 


2007; McCarthy et al. 2010). 


m each of the two objects, we computed the loca¬ 
tion of the spine along the major axis of the filament 
and computed the profiles transverse to it. The spine 
is computed by first identifying the major axis of each 
structure (ba sed on the bounding box of each object, 
as in Section [XT] ) and then computing the centre of gas 
mass in consecutive slices (each one cell thick) perpen¬ 
dicular to the main axis, as shown in Figure [16] Fig¬ 
ures 1 1 7H 1 8| show the profiles of gas density, tempera¬ 


ture ana enclosed gas mass transverse to the main axis 
of the filaments, out to rgi ~ 3.3 Mpc (FI) and 2 Mpc 
(F2). Since the transverse size of each object might 
vary along the spine (Figure |T6|) our radial bins are nor¬ 
malised to the transverse size of the filament along the 
spine, i.e. we bin our profiles according to r/rgi, where 
r-fi | is defined as the minimum distance of the filament 
edge from the spine in the corresponding transversal 
plane. The density profiles for the same physical model 
are rather similar. In the non-radiative case, the cen¬ 
tral mass density is several ~ 1CU 28 g/cm 3 becoming 
~ 2 — 3 times larger in the radiative case, due to the 
higher compression. The mass density smoothly drops 
by a factor ~ 5 — 10 at 0.5rgi from the spine of each fil¬ 
ament, then slowly declining outwards. Notice that the 
density at the last radial bin is larger than our fiducial 
threshold density (agi) due to the averaging at r = rgi 
which can include cells with mass density above the 
threshold resulting from the compression produced by 
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FI &F2 



FI & F2 



Figure 17. Transverse profile of gas density (left) and gas temperature (right) for the filaments FI (solid) and F2 (dashed), in the non-radiative run and in the 
cooling and feedback run (c2). The horizontal lines in the first panel show the density threshold corresponding to different values of \ (Section [33}. The 
radial coordinate is normalised to the transversal size of each filament along the spine (see Section |43) . 


strong external shocks. In all cases the central temper¬ 
ature is around ~ 10 6 K. However, the corresponding 
radial profiles are remarkably different. The FI non- 
radiative filament shows a gently increasing tempera¬ 
ture up to its external boundaries where most of shock 
heating takes place. The same FI object is highly af¬ 
fected by AGN feedback and cooling, which leads to a 
higher central over-density and to a higher temperature 
peak ~ O.o/'fii from the spine of the filament. The effect 
is less evident in filament F2. The horizontal lines show 
the mass cut-off induced by different settings of the agi 
parameter. It is clear how setting am > 2 significantly 
limits the fraction of gas mass assigned to filaments. 
This is confirmed by Figure [18] that shows the cumu¬ 
lative gas mass for the various filaments and models, 
showing that most of the mass is contained at large dis¬ 
tances from the spine, ^ 0.3 —0.5rgi (i.e. ^ 700 — 1500 
kpc) and how only using am ^ 2 this mass can be fully 
captured by our algorithm. Furthermore, as evident in 
the temperature profiles, we thus trace the filament out 
to the position of the main shock fronts including most 
of the cell that have undergone compression and heat¬ 
ing. 


In the panels of Figure 19 we present the profiles 
of the gas pressure, the CR pressure, the pressure 
ratio and the shock Mach number for the same objects 
and runs. In both objects, the thermal gas pressure is 
dominant at all radii compared to the CR pressure. 
However, the CR pressure becomes ~ 10 — 20% of 
the thermal pressure at large radii in the non-radiative 
runs, and ~ 40 — 80 % of the thermal pressure in the 
cooling and feedback runs. The latter is a combined 
effect of the decreased gas temperature and of the 
stronger shocks (as an effect of additional feedback 
bursts) in most of the filament volume, which yield a 
more efficient injection of CRs from DSA. However, 


inside most of the filament the gas motions are super¬ 
sonic (|Ryu et al. 2008b] |Vazza et al.||2014a|) and that, 
therefore, both gas and CR pressure are subdominant 
compared to the ram pressure of the gas. 

Future radio telescopes (such as SKA or, in the nearest 
future, LOFAR) will be able to investigate this regime 
and might provide new clues about the degree of 
plasma collisionality and cosmic rays in this environ¬ 
ment, by directly imaging the morphology of strong 
accretion shocks in radio. One significant emission 
channel for CRs is represented by direct synchrotron 
radiation by accelerated CR-electrons at strong shocks 
(Brown [2011] |Araya-Melo et al.||2012]|. In this case, 
also the level of large-scale magnetisation will affect 
the chance of detection, and present uncertainties 
about the WHIM leave room for several scenarios 


(Ryu et 

pofe ). 


5 DISCUSSION 

Compared to previous work (see Section [T]) that stud¬ 
ied the distribution of DM in the cosmic web (mostly 
using DM-only simulations), our work focused on 
the thermodynamic characterisation of the baryon 
component of the cosmic web. 

The methods that we have used are particularly sen¬ 
sitive to gas densities close to the average cosmic 
value and have no artefacts caused by a particle-based 
DM distribution that especially affect poorly sampled 
low-density regions. As a consequence, the average 
over-density of our population of filaments is smaller 
than what has been quoted in the literature (e.g. |Col- 
berg et al. |2005 ). 
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Figure 19. Transverse profiles of gas pressure (top left), CR pressure (top right), pressure ratio (bottom left) and shock Mach number (bottom right) for 
filaments FI and F2, for the non-radiative run and in the cooling and feedback runs. 


This procedure enables us to extract the cosmic web 
down to an over-density very close to the critical den¬ 
sity and to segment it into its filamentary components. 
Thus, we can analyse a large dynamic range in the 
mass/volume distribution of filaments, using the recent 
2048 3 and 1 024 3 grids produced by our group (Vazza 
et al.|2014bj). This led to the extraction of a very large 
number of filaments at high resolution, e.g. ~ 30000 
objects in our 300 3 Mpc 3 volume, significantly higher 
than what is typically found in DM-only simulations. 

Our method is simpler th an other algorithms 


Hahn et al. 

2007 Aragon-Calvo et al. 

2007 

Sousbie 

et al. 2008 

Cautun et al. 

2014 Aragon-Ca 

Ivo et al. 


to possible misidentification in several pathological 
cases as, for instance, for irregular filaments (see Sec¬ 
tion 4.1) that in some cases may actually be parts of 
sheets with density comparable to that of filaments. It 
can also be affected by resolution effects, e.g. in small 
clumps whose typical size is close to the resolution 
limit. However, a recent comparison has shown that 
the volume occupied by filaments (and, to a lesser 
extent, the mass locked into filaments) can in principle 
be measured equally well by both density-based and 
more sophisticated methods (Cautun et ah] |2014j). 
We conclude that, while additional checks might be 
necessary to assess the nature of specific low mass/size 
objects (e.g. also accessing the DM properties of small 
objects), in general our method gives robust results. 
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We find that the total volume occupied by fila¬ 
ments ranges from 3.8 to 4.5% of the computational 
volu me, which is c onsiste nt with the range provided 
by |Cautun et al.| (|2014J>, and significantly smaller 
than what is found b y previous N-body simula tions 


(e.g. |Hahn et al.12007^] Aragon-Calvo et al.||20T0) The 
largest objects identified in our dataset have a mass 
and estimated length of the order of what is found by 


large N-body sim ulations (e.g. Colberg et al. 2005 
Cautun et alj|2014l, provided that we rescale our gas 


masses by the cosmic baryon fraction in order to get 
the total mass (which is legitimate here since at the 
typical filaments’ over-density the baryon fraction is 
close to the cosmic baryon fraction). However, the 
mass distribution of filaments in the cosmic volume 
(Figure [T2| ) shows a more extended power-law be¬ 
haviour aown t o the smallest masse s compared to the 
literature (e.g. |Cautun et al.||2014| Figure 53). Here 
we can measure an unbroken power-law distribution 
in the range 10 11 — 10 14 Mq, suggesting that our 
sample of filaments is complete down to total masses 
of 6 • 10 11 Mq, i.e. 1-2 orders of magnitude deeper than 
what is typically achieved with N-body simulations. 


We also performed an analysis of the properties 
of individual filaments, focusing on two massive ob¬ 
jects as representative cases. The profiles of gas den¬ 
sity across the filaments is similar to what is reported 


elsewhere (Colberg 

et al. 

20051 

Aragon-Calvo et al. 

2010) The 


|Dolag et al. 

2006 

description ol 

sharp 


tion shocks is a feature in which our simulations are 
particularly effective. Qualitatively similar results for 
the density and temper ature profiles were reported by 
(Klar & Miicket 2012), but these were based on ide¬ 
alised simulations of a forming filament that were not 
set in a cosmological framework. The pattern of strong 
accretion shocks surrounding filaments are consistent 


with previous works in the literatu re (Ryu et al. 2003; 
Skillman et al. 2008[|Vazza et al.||2009 ' ZOn^ Araya- 


Melo et ai.||2012|), which also predict high CR accel- 
eration efficiency at the scale of filaments (Pfrommcr 


et al. 2006). Based on our run-time modelling of (JR 
physics (Vazza et al. [2014b l, the pressure ratios within 


filaments show that CRs can contribute up to ~ 10 — 20 
percent of the total pressure. However, the present un¬ 
certainties CR-acceleration effici encies, especially in 
(low-density) environments (e.g. |Kang & Ryu|[2013j) 
force us to treat this number with caution. The impact 
of radiative cooling and feedback from AGN on the 
WHI M has been explored by several authors (e.g Cen 


& Qstrik e7|2006t |Roncarelli et al.|2006[ [Hallman et al. 


2007). Still, as far as we know, our work represents 

tHenrst study of how cosmic filaments change with 
feedback. However, our modelling of radiative cool¬ 
ing over-estimates the gas metallicity at high redshift, 
and possibly also in the outer parts of filaments at low 
redshifts, even if the observational constraints are still 
weak (e.g. |Urban et al.|2014| ). 


6 CONCLUSIONS 

In this work we investigated the main properties of 
the baryonic matter in cosmic filaments by using a 
large set of cosmological grid simulations recentl y ob¬ 
tained with the ENZO code (Vazza et alj |2014bj ). We 
built a filament identification procedure upon the Visit 
data analysis and visualisation software, exploiting a 
combination of its Isovohtme and Connected Com¬ 
ponents algorithms. The method separates over-dense 
from under-dense regions, identifying a filament as a 
connected set of cells with mass density above a given 
threshold, ogi. The filament identification is further re¬ 
fined by eliminating clusters and small clumps and ac¬ 
cepting only elongated objects. The resulting method¬ 
ology depends on six parameters, the most important 
being the mass density threshold am, whose value lies 
in the range 0.5 — 2.0. The exact value can influ¬ 
ence some statistical properties of the identified objects 
(e.g. the maximum size of filaments extracted from a 
simulation). Other properties (e.g. the average density 
within filaments) are unaffected. The remaining param¬ 
eters are either more tightly constrained or they have 
only a minor impact on the results (see Appendix A for 
further discussion). 

Our most significant findings can be summarised 
as follows: 

• Morphology and thermodynamic properties: fila¬ 
ments, especially long objects (longer than ~ 7 Mpc), 
show a broad variety of shapes and thermodynamic 
properties, that depend on the environment in which 
they lie and on the evolutionary stage. A first rough vi¬ 
sual classification as been attempted, but further inves¬ 
tigation is needed to reach comprehensive and robust 
conclusions. 

• Average properties of the WHIM: the population 
of filament extracted setting q ^ qq and removing col¬ 
lapsed halos, is of the order of ~ 3500 objects within 
(100 Mpc) 3 . The enclosed gas density averaged over 
the whole population is ~ 3 — 5 po, the average tem¬ 
perature is a few r\_/ 10 5 K. 

• Scaling relations: filaments follow well-defined 
scaling relations in the (T, M) and ( V., Al) plane, with 
slopes similar to those of galaxy clusters but differ¬ 
ent normalisation (i.e. lower temperature and larger 
masses). The observed scaling shows that also in fil¬ 
aments gravity sets a clear dependence between the 
enclosed mass and the gas temperature, (T oc M 2 / 3 ) 
even if at low masses the effect of radiative cooling can 
steepen the relation significantly, similar to the case of 
radiative galaxy clusters. 

• Massive filaments: the most massive objects 
found in our suite of simulations have gas masses in 
excess of ~ 1O 15 M 0 , an average (volume-weighted) 
temperatures of ~ 10‘ K and a total volume of a few 
~ 10 3 Mpc 3 , reaching a total length of the order of 
100 Mpc. Only about 1 object of such size can be 
found on average within 100 3 Mpc 3 . 

• Smallest filaments: the smallest filaments that our 
algorithm can reliably extract have a length of about 1- 
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2 Mpc, gas mass of a few ~ 10 10 M© and temperature 
of ~ 10 4 — f 0 5 K (or down to ~ 10 3 K in the radiative 
case). 

• Resolution effects: spatial resolution can affect 
the population of filaments at all scales. The most evi¬ 
dent effect is the high diffusion in low resolution sim¬ 
ulations, that leads to the percolation of otherwise dis¬ 
connected structures. This leads to a drop of the num¬ 
ber of small objects per unit of volume and to an in¬ 
crease in the number of large objects. As a conse¬ 
quence, the same trend can be found in the mass distri¬ 
butions. Temperature distributions show a similar be¬ 
haviour due also to the largest distance travelled by 
shock waves in low resolution simulations, heating up 
larger volumes. In summary, our analysis on the ther¬ 
modynamical and statistical properties of the WHIM 
suggest that a spatial resolution equal or better than 
~ 100 kpc/h is necessary to have converging results 
on the simulated filaments properties, considering ob¬ 
jects of masses larger than 10 10 Mq. 

• Effects of radiative cooling: cooling (mostly dom¬ 
inated by line cooling here, given the assumed fixed 
large metallicity) decreases the average temperature of 
the WHIM by a factor ~ 3 — 5 compared to the non- 
radiative case. In the presence of small-scale clumps 
contained in filaments, the minimum temperature in¬ 
side filaments is reduced by more than one order of 
magnitude. The gas density profile is also increased, by 
a factor ~ 3 close the axis of filaments in the presence 
of cooling. 

• Effects of AGN feedback: efficient feedback from 
AGN can affect most of the volume of filaments, even 
out to large distances from the nearby halos. The ther¬ 
mal feedback implemented in our model causes the 
expulsion of gas from halos and the expansion of fil¬ 
aments, that can percolate and produce significantly 
larger objects compared to the non-radiative case. 

• Effects of cosmic rays: the dynamical impact of 
CRs is small since where the gas density is higher, 
close to the axis of filaments, the ratio between CRs 
and the thermal gas is small (~ 10 — 20 percent in the 
non-radiative case). However, the combination of cool¬ 
ing and AGN feedback is found to increase the budget 
of CRs significantly outwards, where strong accretion 
shocks dominate. 


As a final remark, we stress that a comprehensive 
description of the observable properties of the cosmic 
web requires further important physical mechanisms 
to be taken into account, such as star formation and 
chemical enrichment of the diffuse medium by galac¬ 
tic activity (which affects the distribution of elements 
in the WHIM, and determines the emission/absorption 


properties of the ga s through lines, Cen & Chisarue.g. 


2011 Nicastro et al. e.g. 20131, partial thermal equi- 


bration (e.g. Rudd & Na 

gai [2009), thermal conduc- 

tion (Klar & Miicket 20 

2) and large-scale magnetic 

fields (e.g. Briiggen et al. 

2005; 

Dolag et al. 2008; [Ryu 

et al. 2008a; Vazza et al. 

2014a 

). Such extensive stud- 

ies will be the subject of forthcoming work. However, 


the methodology presented in this paper does not de¬ 


pend on the physics of the problem and can be applied 
to many kind cosmological simulations as is based on 
simple assumptions. It adapts very well to even larger 
datasets, and its parallel strategy will allow us to ex¬ 
ploit it effectively for the identification and character¬ 
isation of filamentary structures in any future simula¬ 
tions programme. 
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APPENDIX A: PARAMETERS AND PERFORMANCE 
ANALYSIS 


Parameters analysis 


Our filaments identification procedure is characterized 
by six parameters, the most relevant being the mass 
density threshold am- Its value has been set accord¬ 
ing to various qualit ative and statistical indications, as 
described in Section 13.31 Although agi = 1 has been 
set, alternative values in the range 0.5 — 2.0 are ac¬ 
ceptable. In order to verify the influence of am on the 
results, we have calculated for simulation 1-1 _1024 
the Mass-Volume relation obtained adopting am = 
0.5,0.75,1.0, 2.0. Figure[20|compares the distributions 
obtained for the different values of am, showing that 
the trends and the dispersions obtained for the differ¬ 
ent thresholds are similar, hence the obtained scaling 
relations are reliable. On the other hand, the mass of 
the largest objects increases with decreasing am, due 
to percolation at low thresholds, and the normalizations 
change. This indicates that the methodology is robust 
for comparative studies of different models or multiple 
realisations of the same model. Absolute estimates (as, 
for instance, the number of objects with a given mass) 
have to be taken with some care. 

The remaining five parameters are: 


(i) a c p mass density threshold to select the higest 
peaks of the mass distribution; 

(ii) 6: scale factor to get clusters’ radii of the order 
of 1 Mpc; 

(iii) Vj. es : objects whose volume is below this value 
are discarded; 

(iv) a: maximum ratio between the axes of the 
bounding box containing a filament; 

(v) ip\ filling factor; filaments with ratio between 
their volume and the associated bounding box volume 
is above this value are discarded. 


The first three parameters can be considered as 
fixed, being related to the physical properties of galaxy 
clusters (the first two) and to the spatial resolution of 
the simulations (the third). For the remaining two pa¬ 
rameters, we have performed the identification proce¬ 
dure several times, changing their values in order to 
investigate the influence on the results. Again, we have 
used simulation 1-1 1024 as a representative case. 


Figure [21] shows the Mass-Volume relation for all 
twelve cases obtained setting a = 1,2,4 and = 
0.1,2.5,5,10, this last parameter being the ratio be¬ 
tween the bounding box diagonal len gth a nd the radius 
of the inscribed cylinder (see Section [33] for more de¬ 
tails): the larger the value of L v , the smaller the thresh¬ 
old filling factor. Figure [22] show a subset of the mass- 
volume data, obtained selecting only those values ob¬ 
tained for Lp = 5 and a = 2 respectively (the val¬ 
ues adopted in this work). Figures [23] and [24] show the 
Mass-Temperature relation for the same parameter set¬ 
tings. 

It is clear how changing the value of these two 
parameters has a minor impact on the results, all the 
main characteristics of the presented statistics being 
preserved. However, a proper selection of their values 
allows eliminating outliers and objects whose geome¬ 
try is not compatible with that expected for typical fil¬ 
aments (e.g. small, round shaped clumps), fine tuning 
the results. The values = 5 and a = 2 have proved 
to be effective without being too restrictive, leading to 
the loss of significant objects. 


Performance analysis 

The data analysis has been performed on the Pilatu^ 
HPC data processing system of ETHZ-CSCS, using 
the client-server 2.8 MPI-parallel version of Visit. Pi- 
latus is a cluster composed of 44 nodes each with 2 
Intel Xeon CPUs E5-2670 (2.60GHz, 16 cores, 64GB 
RAM). 

In a first series of tests, strong scalability (i.e. the 
change in computing time with increasing the number 
of MPI tasks) is investigated. The 2-11024 model is 
adopted as reference. Table [4] presents the CPU time 
as a function of the number of MPI tasks. Accept¬ 
able scalability (though not linear) is achieved up to 16 
tasks, although at 32 and 64 tasks the computing time is 
still decreasing. In the same table, the memory usage is 
shown. Due to the data replica needed for parallel pro¬ 
cessing, the memory required grows with the number 
of tasks, reaching, at 64 tasks, more than three times 
that required for a single task, which is about 17GB. 
Two variables, mass density and temperature, are used 
each represented by a 1024 3 cells floating point single 
precision variable, accounting for 8 GB. The remain¬ 
ing 50% of memory is required to handle isosurfaces 
and connected components. 

For the entire data analysis, we have set the number 
of MPI tasks to 32. This represents a good trade-off 
between computing time and memory usage. In Table[5] 
we present the CPU time needed to run the filament 
reconstruction procedure for the different simulations. 
Only the 1-1_2048 case required 64 MPI tasks, due to 
memory constraints. The “0”, “cl” and “c2” models 
have similar requirements. 

The usage of 32 tasks keeps the CPU time reason¬ 
ably low. This is particularly important during the set¬ 
up of the different model parameters when the filament 


° http://user.cscs.ch/computing_resources/pilatus 
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Figure 20. Mass-Volume relation for simulation 1-1.1024 for different values of a^\. 
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Figure 21. Mass-Volume relation for simulation 1-1.1024 for a = 1.0 (red), a = 2.0 (blue), a = 4.0 (green) and L^ = 0.1 (crosses), L^ = 2.5 
(diamonds), L^ = 5.0 (triangles), L^ = 10.0 (squares). The overall influence of the different settings on the Mass-Volume statistics is shown. 


identification procedure must be repeated many times, 
with slightly different values of each parameter. Note 
that the analysis of the I -1 2048 case requires a huge 
memory allocation (more than 185 GB). Furthermore, 
in the 1024 3 models, the memory usage grows from 3- 
1 _512 to 1-1 _2048 runs because of the increasing simu¬ 
lation volume and hence the larger number of identified 
objects. 


Finally, we analysed the load balancing among 
MPI tasks, to check if the workload is evenly dis¬ 
tributed among the cores. Again, we adopt as reference 
the 2-l_1024 model. As a measure of the load balanc¬ 
ing, in Table 16] we show the minimum and maximum 
memory useaby the various tasks. This can be as¬ 
sumed as an effective measure for the workload, due to 
the homogeneous spatial distribution of the filaments. 
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Figure 22. Mass-Volume relation for simulation 1 -1 _1024 for a = 1.0 (red), a = 2.0 (blue), a = 4.0 (green) and = 5.0 (left panel) and L ^ = 0.1 
(black), Lcp = 2.5 (red), L^ = 5.0 (blue), L ^ = 10.0 (green) and a = 2.0 (right panel). The impact of the various settings for each of the two parameters 
individually on the Mass-Volume statistics is shown. 



Mass 


Figure 23. Mass-Temperature relation for simulation 1-1.1024 for a = 1.0 (red), a = 2.0 (blue), a = 4.0 (green) and L= 0.1 (crosses), L^ = 2.5 
(diamonds), L^ = 5.0 (triangles), L^ = 10.0 (squares). The overall influence of the different settings on the Mass-Temperature statistics is shown. 


Table 4. Strong scalability for the 2-1.1024 model. Computing time (col¬ 
umn 2) and memory usage (column 3) are shown as a function of the num¬ 
ber of cores (column 1). 


MPI tasks 

CPU time (sec) 

Memory (GB) 

1 

633 

17.309 

2 

369 

18.381 

4 

207 

19.296 

8 

123 

21.229 

16 

83 

25.256 

32 

63 

34.592 

64 

50 

53.158 


Table 5. Computing time (column 2) and memory usage (column 3) for 
selected models using 32 cores.) 


Model 

CPU time (sec) 

Memory (GB) 

3-U512 

20 

18.166 

3-1.1024 

52 

32.526 

2-1-1024 

63 

34.592 

1-1.1024 

67 

36.785 

1-172048 

771 

185.541 


In all tests, the difference between the maximum and 
minimum values, so the work imbalance, is between 
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Figure 24. Mass-Temperature relation for simulation 1-1 _1024 for a = 1.0 (red), a = 2.0 (blue), a = 4.0 (green) and L^ = 5.0 (left panel) and L^ =0.1 
(black), Lcp = 2.5 (red), L^ = 5.0 (blue), L^ = 10.0 (green) and a = 2.0 (right panel). The impact of the various settings for each of the two parameters 
individually on the Mass-Temperature statistics is shown. 



Table 6. Minimum (column 2) and maximum (column 3) memory usage on 
the different MPI tasks for the 2-1.1024 model as a function of the number 
of cores (column 1) 


MPI tasks 

Min Mem (GB) 

Max Mem (GB) 

4 

4.553 

5.140 

8 

2.586 

2.817 

16 

1.443 

1.703 

32 

1.011 

1.151 

64 

0.770 

0.868 


the 11 and 15 %. 

The results of the performance analysis clearly show 
the feasibility of this approach for big datasets by ex¬ 
ploiting HPC systems and parallel processing, and its 
potential scalability to the analysis of future larger sim¬ 
ulations. 
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